/*
 * Copyright 2011 Marcin Blazejewski
 *
 * This file is part of Ground0.
 *
 * Ground0 is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * Ground0 is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with Ground0. If not, see <http://www.gnu.org/licenses/>.
 *
 */

/**
 * @file   geometry.hpp
 * @brief  Geometry storage classes and intersection detection 
 * algorithms.
 * 
 */


#ifndef GEOMETRY_H_
#define GEOMETRY_H_

#include <list>
#include <vector>
#include <string>

using namespace std;

/** 3D vector class 
 */
class Vector {
public:
    //@{ axis vectors
    static Vector x_axis;
    static Vector y_axis;
    static Vector z_axis;
    //@}

    //@{ coordinates defining the vector
    float x, y, z;
    //@}

    /// create a new instance with given coordinates
    Vector(float _x=0, float _y=0, float _z=0);

    //@{ standard operators   
    Vector& operator +=(const Vector &v);
    Vector& operator -=(const Vector &v);
    Vector& operator *=(const float f);
    Vector& operator /=(const float f);
    Vector  operator +(const Vector &v);
    Vector  operator -(const Vector &v);
    Vector  operator *(const float f);
    Vector  operator /(const float f);
    bool operator ==(const Vector &v);
    //@}

    /// Rotate the vector by a unit quaternion
    Vector& operator *=(const class Quaternion &q);
    Vector  operator *(const class Quaternion &q) const;

    /// Multiply 3x3 matrix by this vector
    Vector  left_multiply(const float m[][3]) const;

    /// dot product of this vector and vector r
    float   dot(const Vector &r);
    /// cross product of this vector and vector r
    Vector  cross(const Vector &r);
    
    /// return the lenght of the vector
    float                   length();
    /// make this vector a normalized vector (length = 1)
    void                    normalize();
    /// return a normal vector from this vector, keep this object unchanged
    Vector                  get_normal() const;
    /// return squared distance to a point
    float                   distance2_to(const Vector &v) const;
    
    /// create a vector from polar coordinates (radians)
    static Vector           from_polar(float r, float head, float elev);
    /// set this vector to values from polar coordinates (radians)
    Vector&                 set_polar(float r, float head, float elev);
};


/** 
 * Calculate normal vector for CCW triangle defined by v1, v2, v3
 * 
 * @param v1 trianglee definition vertex
 * @param v2 trianglee definition vertex
 * @param v3 trianglee definition vertex
 * @param n return: face normal vector
 * 
 * @return n
 */
Vector& face_normal(Vector &v1, Vector &v2, Vector &v3, Vector &n);

/// legacy wrapper for reference-return version
Vector face_normal(Vector &v1, Vector &v2, Vector &v3);

/** 
 * Multiply matrix m by vector v
 */
Vector operator *(const float m[16], const Vector &v);

std::ostream& operator<<(std::ostream&, const Vector &v);


/** Triangle face class.
 */
class Face {
public:
    /// vertices defininig the triangle
    Vector                  *vertices[3]; 
    /// normals for each vertex
    Vector                  *normals[3]; 
    /// texture coordinates for each vertex
    Vector                  *texcoords[3]; 

    /// Create new face with vertex indecies v0, v1, v2
    Face(Vector *v0, Vector *v1, Vector *v2);
    /// copy contructor
    Face(const Face&);
    /// assignment operator
    Face&                   operator =(const Face&);
};

ostream& operator<<(ostream&, const Face &f);

/******************************************************/
/* MATRIX */
/******************************************************/

float mat_det_3x3(const float A[][3]);
void mat_inverse_3x3(const float A[][3], float ret[][3]);

ostream& operator <<(ostream&, const float A[][3]);

/******************************************************/
/* ALGORITHMS */
/******************************************************/

/** Distance of point P to plane (V0, Vn)
 * Distance of point P to the plane is calculated as the dot product
 * of vector from V0 to P with the plane normal vector.
 * @param the point in space
 * @param V0 any point on plane
 * @param Vn plane normal vector, must be normalized
 * @returns distance of point P to plane defined by V0, Vn
 */
float dist_point_plane(Vector &P, Vector &V0, Vector &Vn);

/** 
 * Check if point a is between points b1,b2.
 * All points must lie on the same line.
 */
bool is_overlap(Vector &a, Vector &b1, Vector &b2);

/** 
 * Return the distance of point P to plane (V0,Vn)
 */
inline float dist_point_plane(Vector *P, Vector *V0, Vector *Vn) {
    return Vn->dot(*P - *V0);
}

/** 
 * Find intersection point between ray P0,P1 and plane V0,Vn
 * Ray is defined as starting at P0 and going through P1.
 * Plane is defined by any point on the plane V0 and normal vector Vn.
 *
 * @param P0 starting point of the ray
 * @param P1 any other point on the ray
 * @param V0 point on the plane
 * @param Vn plane normal vector (must be normalized)
 * 
 * @return intersection point
 */
Vector intersect_ray_plane(Vector *P0, Vector *P1, Vector *V0, Vector *Vn);

/** 
 * Find segment that is the itersection between triangle (tv1,tv2,tv3) and plane (V0,n)
 * 
 * @param tv1 triangle vertex
 * @param tv2 triangle vertex
 * @param tv3 triangle vertex
 * @param V0 plane point
 * @param Vn plane normal
 * @param i1 return segment point 1
 * @param i2 return segment point 2
 * 
 * @return true, if given objects intersect
 */
bool intersect_triangle_plane(Vector *tv1, Vector *tv2, Vector *tv3, 
                              Vector *V0, Vector *Vn,
                              Vector &i1, Vector &i2);

/** 
 * Find intersection segment between two triangles
 * 
 * @param tv11 first triangle, vertex 1
 * @param tv12 first triangle, vertex 2
 * @param tv13 first triangle, vertex 3
 * @param tn1 first triangle normal vector (must be normalized)
 * @param tv21 second triangle, vertex 1
 * @param tv22 second triangle, vertex 2
 * @param tv23 second triangle, vertex 3
 * @param tn2 second triangle normal vector (must be normalized)
 * 
 * @return true, if the triangles intersect
 */
bool intersect_triangle_triangle(Vector *tv11, Vector *tv12, Vector *tv13, Vector *tn1,
                                 Vector *tv21, Vector *tv22, Vector *tv23, Vector *tn2);


/** 
 * Find intersection segment between two triangles
 * 
 * @param tv11 first triangle, vertex 1
 * @param tv12 first triangle, vertex 2
 * @param tv13 first triangle, vertex 3
 * @param tn1 first triangle normal vector (must be normalized)
 * @param tv21 second triangle, vertex 1
 * @param tv22 second triangle, vertex 2
 * @param tv23 second triangle, vertex 3
 * @param tn2 second triangle normal vector (must be normalized)
 * @param i return: intersection segment
 * 
 * @return true, if the triangles intersect
 */

bool intersect_triangle_triangle(Vector &tv11, Vector &tv12, Vector &tv13, Vector &tn1,
                                 Vector &tv21, Vector &tv22, Vector &tv23, Vector &tn2,
                                 Vector i[2]);

// debug version of the above
bool d_intersect_triangle_triangle(Vector *tv11, Vector *tv12, Vector *tv13, Vector *tn1,
                                   Vector *tv21, Vector *tv22, Vector *tv23, Vector *tn2,
                                   Vector i[2]);

#endif
